Extracellular DNAses Facilitate Antagonism and Coexistence in Bacterial Competitor-Sensing Interference Competition

ABSTRACT Over the last 4 decades, the rate of discovery of novel antibiotics has decreased drastically, ending the era of fortuitous antibiotic discovery. A better understanding of the biology of bacteriogenic toxins potentially helps to prospect for new antibiotics. To initiate this line of research, we quantified antagonists from two different sites at two different depths of soil and found the relative number of antagonists to correlate with the bacterial load and carbon-to-nitrogen (C/N) ratio of the soil. Consecutive studies show the importance of antagonist interactions between soil isolates and the lack of a predicted role for nutrient availability and, therefore, support an in situ role in offense for the production of toxins in environments of high bacterial loads. In addition, the production of extracellular DNAses (exDNases) and the ability to antagonize correlate strongly. Using an in domum-developed probabilistic cellular automaton model, we studied the consequences of exDNase production for both coexistence and diversity within a dynamic equilibrium. Our model demonstrates that exDNase-producing isolates involved in amensal interactions act to stabilize a community, leading to increased coexistence within a competitor-sensing interference competition environment. Our results signify that the environmental and biological cues that control natural-product formation are important for understanding antagonism and community dynamics, structure, and function, permitting the development of directed searches and the use of these insights for drug discovery. IMPORTANCE Ever since the first observation of antagonism by microorganisms by Ernest Duchesne (E. Duchesne, Contribution à l’étude de la concurrence vitale chez les microorganisms. Antagonism entre les moisissures et les microbes, These pour obtenir le grade de docteur en medicine, Lyon, France, 1897), many scientists successfully identified and applied bacteriogenic bioactive compounds from soils to cure infection. Unfortunately, overuse of antibiotics and the emergence of clinical antibiotic resistance, combined with a lack of discovery, have hampered our ability to combat infections. A deeper understanding of the biology of toxins and the cues leading to their production may elevate the success rate of the much-needed discovery of novel antibiotics. We initiated this line of research and discovered that bacterial reciprocal antagonism is associated with exDNase production in isolates from environments with high bacterial loads, while diversity may increase in environments of lower bacterial loads.

In order to determine the relative number of antagonists in a soil population, we first tested our ability to retrieve expected ratios of antagonists in a mock community made of E. coli MC4100 (1) and Pseudomonas CVAP#3 (2). To this end, o/n cultures of E. coli and Pseudomonas CVAP#3 were grown at 25˚C in a rotator at max speed (Innova 4130). At the onset of the experiment, E. coli and Pseudomonas CVAP#3 were mixed 1:1 and diluted till a desired dilution in PBS, after with 10-fold dilution series were plated on 10% TSA. After 36h at 25˚C, colonies were counted and their ability to inhibit Staphylococcus CWZ226 (2) determined. As the data in Figure S2-1 indicates, we were able to retrieve the expected number of antagonist from these mockpopulations since no correlation with cell-density was seen in the 1:1 co-inoculated populations. In addition, E. coli and Pseudomonas CVAP#3 were also mixed in different ratios after which 10-fold dilution series in PBS were plated on 10% TSA. After 36h at 25˚C, colonies were counted and their ability to inhibit Staphylococcus CWZ226 determined. Strict correlation (R 2 =0.999) was observed with the different ratios of the mixed E.coli and Pseudomonas CVAP#3. These results support the validity of the method to determine the relative number of antagonists in a population.

Supplemental Material #4A: Development of a Probabilistic Cellular Automaton Introduction
A Cellular Automaton (CA) consists of a toroidal regular planar grid of virtual-cells which evolves through discrete time steps and discrete states, according to rules based on the state of neighboring cells. If the rule is probabilistic, as in our case, the model is a Probabilistic Cellular Automaton (ProbCA) in which a randomly chosen neighbor cell will replace the focal cell according to conditions defined in an interaction matrix. A (Probabilistic) Cellular Automaton is a model useful for simulating complex and nonlinear dynamics and is used in studies to understand bacterial communities. For example, in-vitro and in-vivo game-theoretical dynamics among three phenotypes of bacteria: killer (aggressor), resistant, and sensitive, lead to a stable dynamic equilibrium in the stationary phase and was modeled using a cellular automaton (1,2,3,4). Cellular automata can also help to study coexistence in a dynamic equilibrium and are thus employed to explore the effect of bacterial interactions on biodiversity (5,6,7). For example, using bacterial populations derived from soil, Abrudan (6) used species observed (Richness) as an index to measure the coexistence among Streptomyces. In another study, using isolates from water, an environment less structured than soil, Zapien-Campos (7) found that more aggressive isolates were more likely to occupy the cells in the dynamic equilibrium.
In order to study the effect of antagonist populations from different soils on coexistence, we first developed a Probabilistic Cellular Automaton in Mathematica (11.3.0.0) to explore its use in modeling antagonist interactions within a population of soil derived-bacterial isolates. Mathematica is used in modelling of bacterial interactions (8). Our cellular automaton is an amended version of the model created by Abrudan (6), and we tested its robustness with our experimental data in relation to the following questions: (i) How does changing the grid size affect the outcome? (ii) Is richness to estimate diversity in the dynamic equilibrium a good index? Finally, (iii) does spatial configurations affect the diversity of soil bacteria?

Results and Discussion
The interaction matrix: As mentioned, our ProbCA model of microbial interactions is inspired by the one used by Abrudan (6). The main differences are that we (i) only consider asocial interactions, and no social interactions, (ii) used a randomly seeded homologous ethno-sphere, (ii) of 350x350 cells in size, (iii) a randomly choose the focal cell as well as the neighboring cells, (iv) 100 iterations and (v) 100 simulations in parallel. Furthermore, we (vi) use two populations of antagonists from two different types of soil of different size. The Forest SubSurface (FSS) population was composed of 25 isolates and the Grassland Surface (GS) population was composed of 53 isolates.
In addition, (vii) based on three experimental results, binary interaction matrices were created using a consensus majority rule based on the three observations (Table  S4-1) and rather than the four observations employed by Abrudan et al. (2015). All isolates from a single source were tested pairwise against all isolates from the same source and the ability to produce a zone of inhibition (ZOI) determined ( Figure S4-2A). Development of a ZOI is an indicator of the production of a toxin by the antagonist inhibiting the indicator (or antagonized). There were a total of 2809 tests for the GS population and 625 tests for the FSS population and these tests were repeated three times. Table S4-1: Determination of the consensus in this study. A "1 or +" means that a zone of inhibition (ZOI) was observed on the antagonized by the antagonist. A "0 or -" means that such a ZOI was not observed.

The Probabilistic Cellular Automaton:
For each population we created an interaction matrix as depicted in Figure S4-2A. This interaction matrix was used as input in the simulations. For the simulation, we seeded a square grid randomly with numbers 1 to 25 representing the FSS isolates in a homogenous ethno-sphere ( Figure S4-2C). Every cell on the grid represents an isolate. Every cell was then chosen, in a random order, as a focal cell. One of the eight neighbors of the focal cell was randomly chosen, and if according to the interaction matrix the isolate represented by the neighboring cell produces a ZOI on that of the focal cell, it replaces the focal cell ( Figure  S4-3B). Otherwise both cells retained their isolate. Pseudocode directing this algorithm in Mathematica is provided in Figure  S4-3. Iteration of this process results in the evolution of a dynamic equilibrium ( Figure S4-2C). Like Zapiens-Campos (7), we observed convergence to an equilibrium in all our simulations, as well as clustering of the dominant isolates. We did not quantify either of these observations.

Effect of grid size on the richness of final distribution
Even though Abrudan et al. (6) and Zapiens-Campos et al (7) tested the effect of (a)social interactions of Streptomyces and bacterial isolates from a superficial sediment on coexistence respectively, these authors did not report the effect of different matrix sizes on the final outcome. To test if the grid size affects richness in the dynamic equilibrium, the isolates were used to seed grid sizes of 20x20, 200x200, 350x350 and 500x500, and simulations were run as described earlier. The number of isolates observed in the dynamic equilibrium were low in both Forest Subsurface (FSS) isolates and Grassland Surface (GS) isolates at grid size 20, and were 11.6 and 10.6 respectively. However, as the grid size increases an increase in richness was observed ( Figure S4-4), which plateaued and reached a maximum richness at higher grid sizes. This equilibrium point for the FSS simulation was 24.1 ± 0.2 and the GS simulation results show that it has not reached maximum richness in a dynamic equilibrium yet. Due to this grid-size dependence, we conclude that richness alone is not a reliable parameter to estimate the diversity in the dynamic equilibrium. The reason beyond the richness increases with an increase in grid size can be explained by (i) when the isolates interact on a smaller grid, there is a higher probability of an isolate being eliminated by an antagonist than in a larger grid. This leads to fewer species present in a smaller grid than in a larger grid; (ii) in a larger grid size, the grid may reach a configuration where the sensitive isolate is surrounded by neutral isolates (resistant) and neutral isolates shields attack from antagonists. Moreover, some authors from in-situ studies of microbial diversity mentioned that richness is not a reliable measurement of diversity due to a downward-biased estimator for total species richness and lack of accuracy in measurement (11). In a small sample size, it is hard to observe all species and the species not observed can be under-sampled, thus the richness is overestimated (12). If this also applies to in-silico simulations remains to be seen, but based on these studies and our empirical results, we concluded not to use richness as a measurement for diversity in the dynamic equilibrium alone and thus must be supported

Effect of Grid Sizes on Distribution of the Resulting Community
Because richness alone is not a good predictor of diversity after employing cellular automata, the distribution of the resulting community as a measure for diversity was explored. The prediction is that despite the differences in richness in the dynamic equilibrium, the obtained community structure expressed in frequency would be constant at dynamic equilibrium across different grid sizes. The final frequency distribution was determined at the end of each replication and the mean of the distribution of 100 replicates was calculated. Spearman rank and the Chi 2 test were employed to compare the difference of the frequencies and the rank distributions. The Spearman-rank test indicates that all grid-sizes for GS (rho>0.95, P<3x10 -16 ) and FSS (rho>0.96, P<3x10 -7 ) give the same distributions based on rank (Table S4-Appendix  1A). The Chi 2 test however, has a better ability to detect differences between the frequencies within the distributions. This resulted in a grid-size dependent frequency distribution for the GS populations (0.01<P<0.07) and FSS population (0.23<P<0.25, Table S4-Appendix 1B). The data in Figure S4-5 illustrates these findings. We conclude that for this type of experiment, the Chi 2 test is a too sensitive statistical approach. Spearman, however, suffices.

The Effect of Spatial Configurations on the Community Profile in the Dynamic Equilibrium
We observed convergence to stationary states of our Cellular Automaton after which the cell's spatial configuration no longer changes. Intuitively, this is because the vulnerable strains are shielded by strains that are resistant to the aggressive strain and nonantagonistic towards the vulnerable strains so that the aggressive strains cannot attack the vulnerable strains (7). As in that work, we checked whether the spatial distribution of the populations had an influence on the equilibrium frequency. We did so by iterating our ProbCA until it reached an equilibrium. Usually these equilibria showed substantially clustered populations of dominant isolates (See Figure S4-6). We then randomly scrambled the equilibrium and reran a number of iterations on it, resulting in a new equilibrium. In the case of the FSS, the new equilibrium showed little difference with the scrambled equilibrium, but the dominant isolates increased slightly in numbers, the weaker ones decreased, indicating that the clustering may protect some of the weaker isolates. This effect was substantially stronger for the GS populations. To determine whether the initial spatial configuration of our populations affect the community profile of isolates, we performed statistics on repeated simulations. The final distribution of the consensus experiment using the 350 by 350 grid was utilized and was reseeded randomly for a second round of simulations as shown in Figure S4  The reseeded final distribution was compared with the final distribution of the first simulation. It was hypothesized that the community profile would change due to changing the position of isolates which prevents them from shielding and clustering. The Spearman rank test suggests that both soil populations had a strong correlation between the two final distributions (rho = 0.999 for the FSS simulation and rho = 0.98 for the GS simulation) and the observed correlation was not due to chance (FSS: P=3.3x10 -7 , GS: P=3.2x10 -7 ). However, there was a decrease in richness in both GS and FSS cases (GS ~4 and FSS ~2). The Chi 2 result showed that for the FSS population the difference was not significant. However, for the GS population there was a significant difference between the distributions (FSS: P=0.237, GS: P= 3.52x10 -5 ). A frequency plot of FSS data ( Figure S4-7) showed that in the FSS population, the community profiles are very similar to each other which was supported in the Chi 2 test and Spearman rank test. However, in the GS population, there was a small frequency difference for some of the isolates. Chi 2 identified this difference. This concluded that the spatial configuration does not affect the overall distribution but can affect the final frequency of each isolate.

Conclusion
In conclusion, we developed a cellular automaton in Mathematica allowing us to estimate coexistence of antagonist soil bacterial isolates. Even though many results were published with richness as a measure for coexistence, in the dynamic equilibrium, diversity measures depicting distributions must be included in the argument since distributions, not richness, are independent of grid-size. Spearman Rank is the more appropriate measure to test for significance. Chi 2 is too sensitive unless a P<0.01 is used as a critical point. Using Spearman Rank Correlation, no significant differences in distributions were observed when one changes grid size or spatial distributions.     This movie is the visualization of a simulation using a 50x50 toroidal regular grid of virtual cells. This grid was randomly seeded resulting in a homologous ethno-sphere. The algorithm randomly chooses the focal cell as well as the neighboring cells and applies the replacement rule as explained in SM#4A. This is followed by 100 iterations and 100 simulations in parallel.

Results
To test for antagonist activity towards Staphylococcus sp. CWZ226 (1) and for DNase activity (2), DNA-plates (Difco cat#263220) were seeded with Staphylococcus sp. CWZ226 as described in the material and methods, and soil isolates were tested for their ability to inhibit Staphylococcus sp. CWZ226 and for exDNase activity. The positive control for exDNase activity is Serratia sp. CWZ222 as determined by 16S sequencing. The positive control for inhibition of Staphylococcus sp. CWZ226 is strain Pseudomonas sp. CVAP#3 (1). The negative control was Staphylococcus sp. CWZ226, which should not produce a ZOI, nor harbor exDNase activity. Plates were incubated for three days at 25˚C and the Zone of Inhibition (ZOI) was determined as well as the exDNase activity. The results depicted in Figure S8